Black soliton in a quasi-one-dimensional trapped fermion-fermion mixture 
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Employing a time-dependent mean-field-hydrodynamic model we study the generation of black 
solitons in a degenerate fermion-fermion mixture in a cigar-shaped geometry using variational and 
numerical solutions. The black soliton is found to be the first stationary vibrational excitation of 
the system and is considered to be a nonlinear continuation of the vibrational excitation of the 
harmonic oscillator state. We illustrate the stationary nature of the black soliton, by studying 
different perturbations on it after its formation. 

PACS numbers: 03.75.Ss, 03.75.Lm 
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I. INTRODUCTION 



The one-dimensional nonlinear Schrodinger (NLS) 
equation in the repulsive or self-defocusing case is usually 
written as [| 



iu t + u xx - \u\ 2 u = 0, 



(1.1) 



where the time (t) and space (x) dependences of the wave 
function u(x,t) are suppressed. This equation sustains 
the following dark and grey solitons 



u[x, t) — rix — ct) exp[— i{<p(x — ct) — fit}], (1-2) 



with 



r 2 {x-ct) = rf-2fsecb 2 [£(x-ct)}, (1.3) 
4)(x~ct) = tan _1 [-2£/ctanli{£(£- ct)}}, (1.4) 
£ = V(2^-c 2 )/2, (1.5) 

where c is the velocity, fx is the chemical potential, and 
?7 is related to intensity. Soliton l|1.2(l having a "notch" 
over a background density is grey in general. It is dark 
if density \u\ 2 = at the minimum. At zero velocity 
the soliton becomes a stationary dark soliton or a black 
soliton: \u(x, t)\ = y/rjtanh[xy/r]~/2]. 

The similarity of the NLS equation Ijl.lD to the mean- 
field Gross-Pitaevskii (GP) equation describing a trapped 
zero-temperature Bose-Einstein condensate (BEC) |3j 
imply the possibility of a stationary dark soliton, or 
a black soliton, in a trapped BEC or in a degenerate 
trapped boson-fermion mixture (DBFM) Q. Dark soli- 
tons represent a local minimum in a trapped BEC and 
have been the topic of both experimental 5] and theoret- 
ical HIHB] investigations. Both DBFMpI and degener- 
ate fermion-fermion mixture (DFFM) y| have been ex- 
perimentally observed and studied theoretically [lfjl 1 1 1| . 
It has been suggested that the black soliton of a harmon- 
ically trapped BEC could be a stationary eigenstate Q 
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of the GP equation 3] as in the case of the trap-less NLS 
equation. We demonstrated that the black soliton in a 
trapped zero-temperature BEC 01 or a DBFM ^| is 
the first stationary vibrational excitation of the nonlin- 
ear dynamical equation for the respective systems. No 
such vibrational excitation of the NLS equation (|l.lfl ex- 
ists in the absence of a harmonic trap. 

In this paper we study the numerical simulation of a 
black or a stationary dark soliton in a trapped DFFM. We 
use a coupled time-dependent mean-field-hydrodynamic 
model for a DFFM and consider the formation of black 
solitons in a quasi-one-dimensional cigar-shaped geom- 
etry using numerical and variational solutions. The 
present model is inspired by the success of a similar model 
used recently in the investigation of collapse [l4|. bright 
and dark solitons in a DBFM and of mixing- 
demixing in a DFFM [lj|- The black soliton of a DFFM 
is demonstrated to be a stationary vibrational excitation 
of the mean-field-hydrodynamic model. 

We suggest a scheme for numerical simulation of a sta- 
tionary dark soliton in trapped DFFM by time evolution 
of the linear harmonic oscillator equation starting with 
the analytic vibrational excitation, while the nonlineari- 
ties corresponding to the mean-field-hydrodynamic equa- 
tions are slowly introduced. The simulation proceeds 
through successive eigenstates of the dynamical model. 
To illustrate the stability of our scheme we study the os- 
cillation of the stationary dark soliton upon application 
of different perturbations [p3| . 

Collective excitations in the form of solitons and vor- 
tices in trapped fermions have also been investigated re- 
cently by Damski et al. ^} an d Karpiuk et al. [l8| . 
However, these authors considered isolated ultra-cold 
fermions and not a realistic DFFM as in the experiments 
and as discussed in this paper. Also they did not demon- 
strate the existence of stable dark solitons with a zero 
at the center of the notch as in the present study. They 
identified grey soliton-like structure with a shallow dip in 
the isolated fermionic density distribution quite distinct 
from the stable fermionic dark solitons in a DFFM as 
noted in this investigation. 

Though the present dark solitons are stable in the 
mean-field formulation, they could be unstable physically 
due to quantum fluctuations flfl ] . The effect of quantum 



2 



fluctuations is lost in the mean-free model and can only 
be studied in a field-theoretic approach. Moreover be- 
ing an excited state they are thermodynamically unsta- 
ble. There have been suggestions about how to excite a 
dark soliton by phase imprinting method [T7L ITsj . The 
dark soliton is the lowest vibrational excitation of the 
mean-field model |12| and there have also been inves- 
tigations about how to attain such excited states [20j |. 
Nevertheless, despite different sug gestions about how to 
excite a dark soliton [rH ITsl l20j |. experiments to date 
have not yet generated a stationary dark soliton. Exper- 
imentally, so far the dark solitons have been unstable |4| . 
However, considering that they are stationary excitations 
of the mean-ncld-hydrodynamic equation their creation 
must be possible at least as a non-stationary dark soliton 
which may turn grey and oscillate before decaying due to 
quantum fluctuations and thermodynamic effects. 

The plan of the paper is as follows. In Sec. II we 
present the mean-field-hydrodynamic model which we 
use in the present study. We also perform a reduction of 
the model equations to a quasi-one-dimensional form in 
a cigar-shaped geometry. In Sec. Ill we present a varia- 
tional analysis of the mean-field equations. In Sec. IV we 
present the numerical results for stationary dark solitons 
and compare with the variational results. The station- 
ary nature of the dark solitons is illustrated through a 
numerical study of their dynamics under different per- 
turbing forces. Finally, we present a summary in Sec. 
V. 



II. MEAN-FIELD-HYDRODYNAMIC MODEL 

To develop a set of practical time-dependent mean- 
field-hydrodynamic equations for a DFFM, we consider 
the following Lagrangian density |13J, [15J 



C 



dt 



6m,- 



5 J J 



(2.1) 



where rrij is the mass of component j(— 1, 2), a com- 
plex probability amplitude, nj = \^j\ 2 the real prob- 
ability density, Nj = j drrij (r) the number of respec- 
tive atoms, Aj — fi 2 (Qn 2 ) 2 ^ /(2m,-), and Vj is the 
confining trap. Here the interspecies coupling is gi2 = 
2nTi 2 ai2/mR with the reduced mass uir — mxmz/irni + 
7TI2), and ai2 is the interspecies scattering length. The in- 
teraction between intra-species fermions in spin-polarized 
state is highly suppressed due to the Pauli blocking term 
3AjH^ 3 /5 and has been neglected in Eq. I|2.1I) and will 
be neglected throughout. The kinetic energy terms in 
this equation % A 2 / (fkrij) are derived from a hydro- 
dynamic equation for the fermions [Tlj and contribute 
little to this problem compared to the dominating Pauli 



blocking term in Eq. (|2.1(l . However, the inclusion of 
the kinetic energy terms in Eq. I|2.1|l leads to a smooth 
solution for the probability density everywhere E3- To 
keep the algebra simple and without losing generality, 
in our calculation we shall take equal fcrmion masses: 
mi = 7712 = to/3. This simulates well the fcrmion mix- 
tures in two hyperfine states of 6 Li and 40 K atoms ob- 
served experimentally Q. The Lagrangian density of 
each fermion component in Eq. (|2.1() is identical to that 
used in Refs. El El 

The mean-field dynamical equations for the system are 
just the usual Euler-Lagrange (EL) equations 



d dC 



dt 



2. 

dt 



3 

V 

h dx 



dC 



k g 



,M1/ 



dC 

dm 



(2.2) 



where Xk,k — 1,2,3, are the three space components. 
With Lagrangian density i|2.1[l the following EL equa- 
tions of motion are derived in a straight-forward fashion 
El El: 



d n 2 vl 



dt 6m a 



+ Vj(r)+Ajn 2/3 + 



ih— - - — - + VAt) + + Qi2n k Vj = 0, 

(2.3) 

This is essentially a time-dependent version of a similar 
time-independent model used recently for fermions |22| . 
For a stationary state Eqs. 112.311 y ield the same result 
as the formulation of Ref . jllL Y22\ . For a svstem with 
large number of fermions both reduce to fll. Y22\ the well- 
known Thomas- Fermi approximation [3Ll23|: rij = [{fj.j — 
Vj)/Aj] 3 / 2 , with [ij the chemical potential. 

We reduce Eq. 1)2. 3fl to a minimal quasi-one- 
dimensional form in a cigar-shaped geometry where the 
confining trap of anisotropy v has the form Vj(r) = 
^3rrijUj 2 (p 2 + v 2 z 2 ), with p the radial vector and z the 
axial vector. For a cigar-shaped geometry v << 1, we 
consider solutions of Eq. I|2.3|) of the type $>j(r,t) = 

\fNj4>j{x,T)%l}f\p)/\fli where x = z/l, r = tvuj/2, 

I = (i/lowl) and tpj (p) is the circularly symmet- 
ric ground-state linear harmonic-oscillator wave function. 
In the v << 1 limit the quasi-one-dimensional func- 
tions 6Ax, t) satisfy the following dimcnsionless equa- 
tions 161: 



_5 2 

dx 2 



,4/3 



^(1,^=0, k = 1,2, (2.4) 



where Afjk 
9(6 7 riV j /iy) 2 / 3 /5. 



\2a 12 (N k /v)/l, and M 3 



jj 



In Eq. I|2.4|) . the normalization 



condition is given by \<f>j (x, r)\ 2 dx = 1. 

In case when N\ = N2 = N, the two equations l|2.4|) 
become identical and reduce to 

d 2 



,4/3 



(x,t) =0, 



(2.5) 
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which permits the solution <f> — <fii — 4>2, with (3 = J\fj : 
and 7 = Afjk ■ 



III. VARIATIONAL ANALYSIS 

A. Symmetric Case (TVi — N2) 

Next we present a variational analysis of Eq. (|2.5|) 
based on the stationary trial wave function [24| 



4> v {x,t) = A(r)xexp 



_^_ + l 6(T)x2+zc(r) 



(3-1) 

where A is the amplitude, R is the width, b the chirp, 
and c the phase. The Lagrangian density for Eq. i|2.5[l is 
the one-term version of Eq. H2.1[) . e.g., 



The use of Eqs. (EIJ, and (S3 leads to the follow- 
ing differential equation for the width R: 



d 2 R 
dr 2 
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H 


1 




2\ 


/2^R 2 


d 


' 2 


+ 7 


1r 


R? 


^ R2 



-AR- 
1 



2 2 / 3 a/3 



#5/3 97r l/3 



2?r 



2R Z 



4a(3 



(3.10) 

1/3 



irR 2 



(3.11) 



The quantity in the square bracket is the effective poten- 
tial of the equation of motion. Small oscillation around a 
stable configuration is possible when there is a minimum 
in this potential. The variational result for width R fol- 
lows by setting the right-hand-side of Eq. 13.1l)fl to zero 
corresponding to a minimum in this effective potential, 
from which the variational profile for the soliton can be 
obtained using Eq. (fSTTf) . 



_ ih 
~ 2 
+ x 2 n 



dt 





d(j> v 


+ 


dx 



2 5 

(3.2) 



which is evaluated with this trial function and the effec- 
tive Lagrangian L = J_ oa £(<p v )dx becomes 

L = ^^(*P(AR)^ + c + 3 -^ + ^LA 2 R 2 
2 V 4 16V2 



_L + h 2 R 2 + -R 2 
2R 2 2 2 



(3.3) 



where a = (54/125 v ^F)(3/5) 1 / 6 r(13/6) « 0.242269 and 
the overhead dot denotes time derivative. The variational 
Lagrange equations pH ] 



d dL _dL 
dt dq dq ' 



(3.4) 



where q stands for c, A, R and 6 can then be written as 
^-A 2 R 3 = constant = 1. (3.5) 



6 , 2 3 7 A 2 20 „A 4/3 4c 

36 = - — - 6b 2 



6- — a(3- 



i? 4 2\/2 " 3 ^i? 2 / 3 i? 2 ' 

2 , 57^1 2 52 A 4 / 3 4c 

56 = -^- 106 -lv¥- 10 -¥^^73-^ 

(3.7) 

R = 2Rb. (3.8) 

The constant in Eq. I|3.5|l is fixed by the normalization 
condition. Eliminating c from Eqs. (|3.6() and (|3.7|l we 

obtain 

4 A 2 8 /I 4 / 3 

^=^-^+7^-4+-^. (3.9) 



B. Asymmetric Case (iVi / JV 2 ) 

The above variational analysis can be extended to the 
asymmetric case. However, the algebra becomes quite in- 
volved if we take a general variational trial wave function 
with chirp and phase parameters. As we are interested 
mostly in the density profiles, we consider the follow- 
ing normalized stationary trial wave function for fermion 
component j of Eq. (|2.4I) 



rrexp 



2R (r) 



j = \,2. (3.12) 



Using essentially the Lagrangian density i|2.1[l in this case 
we obtain the following effective Lagrangian 



L 



E>2 jy2 



2 



/a2 2 / 3 Njj 3(i? 4 + l) 



{R\ + R 2 f/ 2 ^W 1/3 R 



2/3 

i 



2R 2 



(3.13) 



The variational Lagrange equations (|3.4|) for i?i and R2 
now become 



i? 3 J 



(i?? + i?i) 7 / 2 



2 2 / 3 SaA^- _ 

97r l/3 „5/3 



(3.14) 



Equations (|3.14l) can be solved for the variational widths 
Rj and consequently the variational profile of the wave 
functions obtained from Eq. I|3.12|l . When N\ = N2, in 
Eq. |3~T4^ Njj = fi, N jk = N k j = 7 an d R x = R 2 = R; 
and in this case it is verified that Eq. (|3.14|) yields the 
same variational widths as from the result obtained in 
the symmetric case in Sec. IIIA given by Eq. (|3.10|l : e. 



4 7 j 2''" apa 

R 3 " + R 2 2~V^ + R 5/3 9ttV3 



7 1 



2 2 / 3 8j3a _ 



= 0. 



(3.15) 
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FIG. 1: (Color online) The solito ns \<j>(x)\ of BEC from the 
numerical solution of Eq. 13.16H (denoted num) vs. x for 
(a) n = 10 and (b) n = 100. The corresponding variational 
solutions obtained from Eq. 13.171 (denoted var) are also 
shown. 



C. Black Soliton in a BEC 

By a variational calculation we establish the stability 
of a black soliton in a trapped BEC satisfying the one- 
dimensional GP equation: 



dx 2 



d 



(x,r)=0. (3.16) 



This equation is identical with Eq. (|2.5|) with f3 = and 
7 = n. Hence the variational analysis of Sec. Ill A applies 
to this case and the differential equation for width R is 
given by Eq. (|3~TU|) with ,0 = 0: 



d 2 R 



4 



AR 



(3.17) 



The variational width in this case is given by setting the 
right-hand-side of this equation to zero. 



IV. NUMERICAL RESULTS 

We solve Eqs. <|2.4|l for bright solitons numerically 
using a time-iteration method based on the Crank- 
Nicholson discretization scheme elaborated in Ref. [25| . 
We discretize the coupled partial differential equations 
(|2.4() using time step 0.0004 and space step 0.02 in the 
domain — 16 < x < 16. The second derivative in x is 
discretized by a three-point finite-difference rule and the 
first derivative in time by a two-point finite-difference 
rule. We start a time evolution of Eqs. (|2.4() set- 
ting the nonlinear terms to zero, and starting with the 
eigenfunction of the linear harmonic oscillator problem: 
<j)(x,r) = a/ (2/ y/n)x exp(— x 2 /2) exp(— 3fr). During the 
time evolution the nonlinear terms are switched on slowly 
and the time evolution continued to obtain the final con- 
verged solutions. In addition to solving the coupled equa- 
tions H2.4fl , we also solved the single equation l|2.5|l in the 
symmetric case: Aq — N2 and Eq. (|3.16|l for bosons. 
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FIG. 2: (Color online) The black soliton |0(a?)| from the nu- 
merical solution of Eq. i12 . ~>l (denoted num) vs. x for a DFFM 
with Aq = 10, N 2 = 10, a 12 /l = 0.05, while Nn = N22 ~ 275, 
and Nij — 60. The corresponding variational solution ob- 
tained from Eq. 13.171 (denoted var) is also shown. 



In our numerical study we take I — 1 /im and consider 
a DFFM consisting of two electronic states of 40 K atoms. 
For a BEC we consider 39 K atoms. This corresponds to 
a radial trap of frequency ui = h/(l 2 m) w 500 Hz. Con- 
sequently, the unit of time is 2/u « 4 ms. For another 
fermionic atom the w value gets changed accordingly for 
/ = 1 /im. 



A. Black Soliton in a BEC 

The results for black solitons in a BEC are obtained 
by solving Eqs. (|3.16(l and H3.17J) and are presented in 
Figs. 1 (a) and (b) for n — 10 and 100, respectively, 
where the soliton profiles are displayed. The variational 
profiles lie close to the numerical results and confirm the 
existence of a black soliton in this case. For smaller values 
of nonlinearity, the solitons have a probability distribu- 
tion close to the ansatz given by i|3.1[) . Consequently, the 
variational result is more accurate. For larger nonlineari- 
ties the actual probability distribution obtained from the 
numerical solution deviates from ansatz l|3.1|) as can be 
seen from Fig. 1 (b). By construction the variational 
ansatz (|3.1|l is stationary. The close agreement between 
the variational and numerical solutions demonstrates the 
stationary nature of these dark solitons. Hence, contrary 
to comments in the literature 0, these black solitons of 
a trapped BEC are stationary. Although, we illustrate 
the present simulation for a black soliton for a repulsive 
BEC, stationary dark soliton also exists for an attractive 
BEC El- 



B. Black Soliton in a DFFM 

The numerical results for black solitons in a DFFM 
are obtained by solving Eqs. 12.511 and 12.4fl for the sym- 
metric (Aq = N2) and asymmetric (i\q 7^ A^) cases, 
respectively. The corresponding variational results are 
obtained from Eqs. (|3.15|l and l|3.14|l . respectively. First 
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FIG. 3: (Color online) The black solitons (a) <f>i(x) and (b) 
4>2 (x) from the numerical solution of Eqs. 12.41 (denoted num) 
for a DFFM with Ni = 10, N 2 = 20, a 12 /l = 0.05, while 
Mi « 275, A/22 « 436, A/12 = 120, and A/21 = 60. The cor- 
responding variational solutions (denoted var) obtained from 
Eqs. 13.141 are also shown. 
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FIG. 4: (Color online) The dark soliton profile for the black 
soliton of a DFFM of Fig. 2 with Ni = N 2 = 10 at times r = 
0, 7r/2, 7r, 37r/2, 27t, 5n/2, and 37r after the soliton has been set 
into oscillation by displacing the trap through a distance of 2 
(im. 
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we present the numerical and variational results for the 
symmetric situation with 012// = 0.05, v = 0.1, and 
Ni = N-2 = 10 in Fig. 2. The nonlinearities here are 
quite large, e.g., Njj w 275 and Njk — 60, j 7^ fc = 1, 2. 
The variational and numerical results in this case com- 
pare favorably as can be seen in Fig. 2. 

Next we present the numerical and variational results 
for the asymmetric case with Ni — 10 and N2 = 20. 
0,12// = 0.05, and v — 0.1 corresponding to even larger 
nonlinearities A/11 « 275, A/22 ~ 436, A/12 = 120, and 
A/21 = 60. The variational and numerical results in this 
case for black solitons for components 1 and 2 are exhib- 
ited in Fig. 3 (a) and (b). 

From Figs. 1—3 we see that although the variational 
result is a good quantitative approximation to the exact 
result for small values of nonlinearity (Figs. 1), it ceases 
to be so for larger values of nonlinearity (Figs. 2 and 
3). This is true, as for large values of nonlinearity the 
harmonic-oscillator ansatz for the soliton profile (|3.1(l is 
not a good approximation to the accurate numerical re- 
sult. In other words, the variational calculation confirms 
a stationary black soliton in all cases. A poor initial 
variational harmonic-oscillator ansatz for the soliton is 
responsible for the imprecise variational results for large 
nonlinearities. 

To demonstrate the robustness of the black solitons, 
we introduced a displacement in the notch of the black 



FIG. 5: (Color online) The simple harmonic motion of the 
minimum of the notch of the black soliton(s) of Fig. 2 after 
the soliton(s) have been set into oscillation by displacing the 
trap through a distance of 2 firci. 



soliton upon formation by displacing the center of the 
harmonic trap by a small distance. The dark soliton(s) 
execute harmonic oscillation for small displacements dur- 
ing a moderate time interval. For the soliton of Fig. 
2 with Ni — N2 = 10 the oscillation is illustrated in 
Fig. 4 by exhibiting successive snapshots of the dark 
soliton at positions of maximum displacement at times 
r = 0, 7r/2, 7T, 37r/2, 27T, 57r/2 and 37r for a trap displace- 
ment of 2 /Ltm at time 0. The soliton remains black at 
all times with the minimum at the notch equaling zero. 
The simple harmonic motion of the position of the notch 
is illustrated in Fig. 5. The angular frequency of oscil- 
lation as obtained from Fig. 5 is the trap frequency 500 
Hz. The numerical simulations presented in Figs. 4 and 
5 demonstrate the stationary nature of the black solitons 
under small oscillations. 

To further demonstrate the stationary nature of the 
black solitons, after their formation we increased sud- 
denly the strength of the harmonic trap by 10%: x 2 — > 
1.1a; 2 . The stability of the dark soliton for ^1=^2 = 10 
under this sudden change in the trap frequency is exhib- 
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FIG. 6: (Color online) The dark soliton profile for the black 
soliton of Fig. 2 with Ni = N 2 = 10 at times t = 0, 40 ms, 
80 ms, 120 ms, and 160 ms after the soliton has been set into 
breathing oscillation by increasing the strength of the trap by 
10%. 



ited in Fig. 6 where we plot the soliton profile at different 
times. The soliton executes reasonably stable breathing 
oscillation. To study the breathing oscillation in detail 
we calculate the root-mean-square (rms) size (s) of the 
soliton at different times from a numerical solution of 
the mean-field equation and the result is plotted in Fig. 
7 (a). We also solved the variational equation (|3.10|l . an 
ordinary differential equation by a fourth-order Runge- 
Kutta rule, to calculate the size of the soliton and the 
time variation of the results is plotted in Fig. 7 (b). The 
amplitude and frequency of the variational oscillation are 
in good qualitative agreement with the precise numeri- 
cal results. However, there is quantitative disagreement 
between the two. After multiplying the variational size 
by 0.89, the results in Fig. 7 (b) are in good agreement 
with the numerical result in Fig. 7 (a). The period of 
breathing oscillation from Fig. 7 (a) is 7.3 ms and that 
from Fig. 7 (b) is 7.6 ms. The fair agreement between 
the variational and numerical results of breathing oscil- 
lation illustrates that the variational description of the 
stationary dark soliton is on the right track and confirms 
the stationary nature of the black soliton. 



numerical and variational methods. The existence of sta- 
tionary dark solitons is predicted by variational study 
and confirmed by detailed numerical analysis. The sta- 
tionary nature of the present solitons is demonstrated nu- 
merically through their (a) sustained oscillation initiated 
by a sudden displacement of the harmonic trap, as well as 
(b) breathing oscillation initiated by a sudden increase of 
the strength of the harmonic trap by 10%: x 2 — > 1.1a; 2 . 
Here we used a set of mean-field equations for the DFFM. 
A proper treatment of a degenerate Fermi gas or DFFM 
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160 



FIG. 7: (Color online) The periodic motion of the rms size (s) 
of the dark soliton after suddenly increasing the strength of 
the trap by 10 %: (a) numerical result, (b) variational result. 



should be done using a fully antisymmetrized many-body 
Slater determinant wave function IT(i| as in the case of 
scattering involving many electrons 26]. However, in 
view of the success of a fermionic mean-field model in 
studies of collapse [l4|. bright 0] and dark 01 soliton 
in a DBFM and of mixing and demixing in a DFFM ^fj , 
we do not believe that the present study on bright soli- 
tons in a DFFM to be so peculiar as to have no general 
validity. 



V. SUMMARY 

In conclusion, we used a coupled mean-field- 
hydrodynamic model for a DFFM to study the formation 
of black solitons in a quasi-one-dimensional geometry by 



Acknowledgments 

The work is supported in part by the CNPq and 
FAPESP of Brazil. 



[1] Y. S. Kivshar and G. P. Agrawal, Optical Solitons - From 
Fibers to Photonic Crystals (Academic, San Diego, 2003). 

[2] P. G. Drazin, R. S. Johnson, Solitons: An Introduction, 
(Cambridge University Press, Cambridge, 1989). 

[3] F. Dalfovo, S. Giorgini, L. P. Pitaevskii, and S. Stringari, 
Rev. Mod. Phys. 71, 463 (1999); V. I. Yukalov, Laser 



Phys. Lett. 1, 435 (2004); V. I. Yukalov and M. D. Gi- 
rardeau, Laser Phys. Lett. 2, 375 (2005). 

[4] S. A. Morgan, R. J. Ballagh, and K. Burnett, Phys. Rev. 
A 55, 4338 (1997); T. Busch and J. R. Anglin, Phys. 
Rev. Lett. 84, 2298 (2000) 

[5] J. Denschlag et al, Science 287, 97 (2000); B. P. Ander- 



7 



son et al, Phys. Rev. Lett. 86, 2926 (2001); S. Burger et 
al, Phys. Rev. Lett. 83, 5198 (1999); S. Stock, B. Bat- 
telier, V. Bretin, Z. Hadzibabic, and J. Dalibard, Laser 
Phys. Lett. 2, 275 (2005). 
[6] R. D'Agosta, B. A. Malomed, and C. Presilla, Phys. Lett. 
A 275, 424 (2000); N. G. Parker, N. P. Proukakis, C. F. 
Barenghi, and C. S. Adams, J. Phys. B 37, S175 (2004); 
N. P. Proukakis, N. G. Parker, D. J. Frantzeskakis, C. S. 
Adams, J. Opt. B 6, S380 (2004); P. J. Y. Louis, E. A. 
Ostrovskaya, Y. S. Kivshar, J. Opt. B 6, S309 (2004); V. 
A. Brazhnyi, V. V. Konotop, Phys. Rev. A 68, 043613 
(2003). 

[7] P. G. Kevrekidis, R. Carretero-Gonzalez, G. Theocharis, 
D. J. Frantzeskakis, and B. A. Malomed, Phys. Rev. A 

68, 035602 (2003). 

[8] F. Schreck et al., Phys. Rev. Lett. 87, 080403 (2001); 
A. G. Truscott et al, Science 291, 2570 (2001); Z. Hadz- 
ibabic et al, Phys. Rev. Lett. 88, 160401 (2002); G. Mod- 
ugno et al, Science 297, 2240 (2002); C. Ospelkaus, S. 
Ospelkaus, K. Sengstock, and K. Bongs, Phys. Rev. Lett. 
96, 020401 (2006). 
[9] B. DeMarco and D. S. Jin, Science 285, 1703 (1999); K. 
M. O'Hara et al, Science 298, 2179 (2002). 

[10] K. Molmer, Phys. Rev. Lett. 80, 1804 (1998); R. Roth, 
Phys. Rev. A 66, 013614 (2002); T. Miyakawa, T. Suzuki, 
and H. Yabu, ibid. 64, 033611 (2001); X.-J. Liu, M. 
Modugno, and H. Hu, ibid. 68, 053605 (2003); X.-J. Liu 
and H. Hu, ibid. 67, 023613 (2003); M. Modugno et al, 
ibid 68, 043626 (2003); D. M. Jezek, M. Barranco, M. 
Guilleumas, R. Mayol, and M. Pi, ibid70, 043630 (2004). 

[11] P. Capuzzi, A. Minguzzi, and M. P. Tosi, Phys. Rev. A 

69, 053615 (2004). 



[12] S. K. Adhikari, J. Low Temp. Phys. (in press), 

cond-mat/0510526 
[13] S K. Adhikari, J. Phys. B 38, 3607 (2005). 
[14] S. K. Adhikari, Phys. Rev. A 70, 043617 (2004). 
[15] S. K. Adhikari, Phys. Rev. A 72, 053608 (2005); S.K. 

Adhikari, Laser Phys. Lett. DOI 10.1002/lapl.200400047. 
[16] S. K. Adhikari, Phys. Rev. A 73, 043619 (2006). 
[17] B. Damski, K. Sacha, and J. Zakrzewski, J. Phys. B 35, 

L153 (2002). 

[18] T. Karpiuk, M. Brewczyk, and K. Rzazewski, J. Phys. B 

35, L315 (2002); E. Witkowska and M. Brewczyk, Phys. 
Rev. A 72, 023606 (2005). 

[19] J. Dziarmaga, Phys. Rev. A 70, 063616 (2004); J. Phys. 
B 36, 1217 (2003). 

[20] V. I. Yukalov, E. P. Yukalova, and V. S. Bagnato, Phys. 
Rev. A 56, 4845 (1997); V. I. Yukalov, K. P. Marzlin, 
and E. P. Yukalova, Phys. Rev. A 69, 023620 (2004). 

[21] H. Goldstein, Classical Mechanics, 2nd Edition, (Addi- 
son Wesley, Reading, 1980). 

[22] Z. Akdeniz, A. Minguzzi, P. Vignolo, and M.P. Tosi, 
Phys. Lett. A 331, 258 (2004); P. Capuzzi, A. Minguzzi, 
and M. P. Tosi, Phys. Rev. A 68, 033605 (2003). 

[23] L. Vichi and S. Stringari, Phys. Rev. A 60, 4734 (1999). 

[24] D. Anderson, Phys. Rev. A 27, 3135 (1983). 

[25] S. K. Adhikari and P. Muruganandam, J. Phys. B 35, 
2831 (2002); P. Muruganandam and S. K. Adhikari, ibid. 

36, 2501 (2003). 

[26] P. K. Biswas and S. K. Adhikari, J. Phys. B 33, 1575 
(2000); 31, L737 (1998); 31, 3147 (1998); 31, L315 
(1998); S. K. Adhikari, Phys. Rev. C 19, 1729 (1979); L. 
Tomio and S. K. Adhikari, Phys. Rev. C 22, 28 (1980). 



